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Abstract 



In this paper a simple (2+1) solid-on-solid model of the epitaxial films growth 

based on random deposition followed by breaking particle-particle lateral bonds 

'^ . and particles surface diffusion is introduced. The influence of the critical num- 

^ \ ber of the particle-particle lateral bonds z and the deposition rate on the surface 

Qj ■ roughness dynamics and possible surface morphology anisotropy is presented. The 

roughness exponent a and the growth exponent (i are (0.863,0.357), (0.215,0.123), 

(0.101,0.0405) and (0.0718,0.0228) for z = 1, 2, 3 and 4, respectively. 

1^^ . Snapshots from simulations of the growth process are included. 

\^ • Key words: anisotropy, computer simulations, molecular beam epitaxy, Monte 
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• • ■ 1 Introduction 

> • 

X 

c^ ■ Theoretical modelling of films epitaxial growth — which is growth of an ori- 

ented single-crystal film of one material upon a single- crystal substrate of 
another [1] and when the main microscopic process is particles deposition fol- 
lowed by their diffusion on the surface — may be grouped to continuum and 
discrete approaches [2]. 
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The continuum ones are based on the stochastic differential Langevin equa- 
tions [3,4,5,6] which group models to the various classes of universality asso- 
ciated with the some set of critical exponents. In discrete approach computa- 
tional methods may be applied, with cellular automata technique (CA) [9,10] 
among others. The latter bases usually on the Arrhenius-like kinetics [11,12] 
and requires knowing many parameters, for instance: 

• material and neighbourhood-dependent activation energies for different el- 
ementary processes, 

• adatom/adatom, substrate-atom/substrate-atom and substrate-atom/adatom 
bonds strengths, 

• vibration factors, 

• substrate temperature, 

• incoming particles flux, etc. 

Usually, the values of these parameters are fitted to reproduce some experi- 
mentally measured particles and/or surface characteristics. 

The title of Ref. [10] characterise the CA technique very well, as for CA, both, 
discrete time and space are necessary. The model must include also the rule 
which tells how the states of lattice cells are subsequently updated. Usually, 
CA models for epitaxial growth simulations [13] base on simple or even toy, 
mechanical rules, e.g. random deposition followed by particles relaxation (see 
Ref. [2] for review). In the relaxation process particles often virtually move to 
the nearest-neighbourhoods (NN) sites to check offered there accommodation 
conditions and then chose the best one. Quite often, the subsequent particle 
arrives to the place of its first contact with surface only when previously 
deposited particle migration process has been completed. Such a situation 
corresponds to very low flux of incoming particles and does not meet real 
condition during molecular beam epitaxy experiments. 

Here we would like to present simple CA which looses both: disadvantages of 
non-physical particles virtual movement to the NN sites and low particles flux 
but still does not need material dependent constants. 



2 Model 



Presented here model, as usually, is extension to random deposition model 
(RDM) with additional particles relaxation. The solid-on-solid approximation 
is applied, so the film may be fully characterised by single-valued function 
h{x,y,t) of the film height in planar coordinates {1 < x < L,l < y < L) in 
time t. 



We start our simulation with perfectly flat substrate. Every rL^ time steps 
new jet of ^dcp-^^ particles arrives. Each time step — between subsequent acts 
of the depositions — particles 'sitting' on the column top may diffuse on the 
surface. The only mobile particles are those which currently have less than z^ 
and Zy created particle-particle lateral bonds (PPLB) in x- and y-direction, 
respectively. For isotropic case only one number z guards the particles mobility. 
Active particles and their movement directions are picked up randomly. The 
particles are not allowed to climb on higher levels, but they are able to jump 
down at the terrace edge. The simulation is carried out until a desired film 
thickness d^^.^ has been deposited. 



3 Results 



We characterise surface morphology with some statistical parameters such as 
film thickness Q and surface width w. The film roughness w is defined as the 
film height h standard deviation 
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from the average film thickness 

e{t)=J2h{x,y,t)/L\ 



x,y 



For anisotropic case {z^ ^ Zy and z^Zy ^ 0) few measures e describing quan- 
titatively surface morphology anisotropy were proposed and investigated. 

For zZxZy = and/or r = the RDM results are well reproduced with Poisson 



distribution of film heights and w{t) = JO{t). The film sample for RDM is 
presented in Fig. 1(a). 

3.1 Isotropic case 



For wide variety of discrete models of surface growth the surface width w is 
expected to follow the dynamical Family- Vicsek [17] scaling law: 
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Fig. 1. Samples of film for (a) RDM [zz^Zy = and/or r 
(b) z = 1 [c) z = 2 (d) z = 4, and anisotropic cases (e) Z; 
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Zy = I (g) Zx = I, Zy = 3 (h) z^ = 2, 
for RDM and isotropic cases and ^max 
the deposition rate was given by ^dcp 
of lattices are presented. 



2 (f) zx = 2, 
:y = 3. Average film height Omax = 2.5 [ML] 
= 0.2 [ML] for anisotropic ones. For all cases 
= 0.01 [ML] and r = 10. 30 x 30-large parts 



where a, j3 and (^ are roughness, growth and dynamic exponent, respectively. 
From Eq. (1) one may read surface roughness dynamics: before reaching the 
characteristic film thickness ^sat the surface roughness grows like w{t) oc t^ 
(see Fig. 3(a) and Tab. 1). Then, the roughness saturates on the level Wsat 
depending on the substrate linear size L: w^at oc L" (see Fig. 3(b) and Tab. 
1). These exponent predict also classes of universality [6,7]. However, deter- 
mination of the exponents a, j3 and (^ solely on the calculation of w{L,t) has 
to be treated with caution [8]. 



For large enough number of relaxations, e.g. r > 5, the characteristic layer- 
by-layer oscillations of roughness w{t) occur and Frank- van de Merwe growth 



Table 1 

Roughness a and growth exponents (3 dependence on the critical PPLB number 
z for 0dep = 0.1 [ML] and r = 1. The results are average over A'run independent 
simulations (from several thousand for L = 5 to a few for L = 100). 

z 1 2 3 4 



a 0.863 0.215 0.1005 0.0718 
(3 0.357 0.123 0.0405 0.0228 



(a) L=1000, z=1, 9dep=0.1 [ML], 6^3,=! [ML] 



(b) L=1000, z=2, ede„=0.1 [ML], e^3,=10 [ML] 




10 20 30 40 50 60 70 80 90 100 
t 

(c) L=1000, z=3, ej„=0.1 [ML], S^^=^0 [ML] 



10 20 30 40 50 60 70 80 90 100 

t 

(d) L=1000, z=4, ed,p=0.1 [ML], 9^3,=10 [ML] 




10 20 30 40 50 60 70 80 90 100 
t 



10 20 30 40 50 60 70 80 90 100 

t 



Fig. 2. Time evolution of the surface roughness w{t) for different number of the 
particle relaxations r and different critical values of PPLB (a) z = 1, (b) z = 2, (c) 
z = 3 and (d) z = 4 (L = 10^ 6ldep = 0.1 [ML]). 



mode is observed (see Figs. 2 and 4). The influence of the critical number of 
PPLB z on the roughness and growth exponents for given 6dep = 0.1 [ML] and 
r = 1 are collected in Tab. 1. 

Surface roughness decreases with r (see Fig. 2) but increases with 6'dep (see 
Fig. 5), as expected. 



The fllm samples for z = 1, 2, 4 are presented in Figs, l(b-d). 
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Fig. 3. Dependence of the roughness dynamics w{t) for (a) z = 1, (h) z = 2, (c) 
z = 3, (d) z = A and (e) the saturation level Wsax. on linear lattice size L for ^dcp = 0.1 
[ML] and r = 1. 
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Fig. 4. Time evolution of the surface roughness w{t) for L = 100, z = 4, ^dep = 0.1 
[ML] and different values of r. 
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Fig. 5. Dependence of the roughness dynamics w{t) on ^dcp for L = 1000, r = 1 
and (a) z = 1 (b) z = 4. 6'dep = 0.01, 0.02, 0.05, 0.1, 0.2, 0.5, 1.0, 2.0 from bottom 
to top. 

Table 2 

Different measures of the surface morphology anisotropy e for L = 1000, ^dcp = 0.01 

[ML], 6lmax = 0.2 [ML] and r = 10. 



1 



Zy 2 



3 
1 



3 
2 



£1 0.424 -0.424 0.457 -0.457 0.013 -0.013 
£2 3.05 0.443 3.32 0.408 1.11 1.04 
£3 2.61 2.61 2.61 2.61 3.56 3.56 



3.2 Anisotropic case 



For submonolayer substrate coverage (e.g. 6'max = 0.2 [ML]), and when Zx ^ 
Zy and ZxZy ^ (see snapshots from simulations presented in Figs, l(e-h)) 
we search for quantitative measure of surface morphology anisotropy. One of 
such characteristic — based on the height-height correlation function — was 
proposed in Ref. [14] and employed for investigation of the surface morphology 
thermal evolution in Ref. [15]. 

Here, for quantitative characterisation of the surface morphology anisotropy 
we propose set of e parameters: e\ = {(p^ — (py)/{(px + (py), £2 = (px/(py, and 
£3 = i/A, where (px and (py are x- and y-side of the minimal rectangle which 
totally covers whole cluster, i is the cluster perimeter and A is the the cluster 
area. The periodic boundary condition modification of the Hoshen-Kopelman 
algorithm [16] for cluster labelling were used. The results are averaged over 
all existing clusters. Note, that for isotropic case (px = 4>y yields £1 = and 
£2 = 1- The values of £3 depend on particular clusters shape and their size, 
i.e. £3 = 4/a for square of side a, £3 = 2/r for circle with radius r and 
£3 = Ay/S/a for equilateral triangle of side a. It seems that £1, for which 
ei{zx,Zy) = —ei{zy, Zx), is much more suitable for quantitative morphology 



characterisation than £2 and/or £3. 



4 Summary 



For all values of the critical PPLB z we observe scaling of the surface roughness 
wiL^ t) according to Family-Vicsek law (1) with characteristic exponents given 
in Tab. 1. 

For fixed set of model control parameters the roughness w decreases both with 
parameter r and z. The latter one may be also seen in a very naive way as 
an equivalent of the absolute substrate temperature, as with increase of z the 
particles become more and more mobile (they stick at the place of deposition 
for z = 0, while for z = 4 their moves become completely forbidden only when 
they have four NN). 

The increase of the parameter r corresponds to effective decrease of the incom- 
ing particles flux per unit substrate area and the time step. For lower particles 
flux we observe smoother surfaces, as in between subsequent acts of deposition 
particles may make more moves. The similar effect may be reproduced with 
direct decreasing of the deposition rate ^dep, and thus surface roughness w 
increases with ^^cp- 

For anisotropic case the dependence of ei parameter may be qualitatively 
compared with results of Monte Carlo simulations from Ref. [12], CA approach 
presented in Refs. [14,15] and — at last but not least — with experimental data 
[12]. So, when one identify z with the substrate temperature one may observe 
initial increase of ei followed by its decrease — i.e. ei{z = 1.5) < ei{z = 
2) > ei{z = 2.5), where z = {zx + Zy)/2. We conclude, that Si has similar 
properties as e defined via height-height correlation function in Ref. [15], and 
additionally is antisymmetric (i.e. ei{zx-,Zy) = —ei{zy,Zx)) and vanishes for 
isotropic case {z^ = Zy and/or z^Zy = 0). 

Presented here asynchronous and probabihstic CA 

• does not need virtual movements of the particles to the NN sites to check if 
they offer better accommodation condition (i.e. higher coordination number 
or lower column height), 

• and do not have many activation energies (for different elementary processes 
like surface diffusion or breaking bonds) and lattice vibration factor, etc. 

The latter may be seen as the model disadvantage, however, one can easily 
imagine how to introduce material constants into growth rules to obtain the 
model which lays somewhere in between CA with simple, mechanical rules 



of surface relaxation and full reversible Arrhenius based kinetic Monte Carlo 
model. 

The program for growth process microscopy is available from our web-page 

[18]. 
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